home *** CD-ROM | disk | FTP | other *** search
/ Space & Astronomy / Space and Astronomy (October 1993).iso / mac / TEXT_ZIP / spacedig / V10_1 / V10_123.ZIP / V10_123
Internet Message Format  |  1991-07-08  |  13KB
  1. Return-path: <ota+space.mail-errors@andrew.cmu.edu>
  2. X-Andrew-Authenticated-as: 7997;andrew.cmu.edu;Ted Anderson
  3. Received: from beak.andrew.cmu.edu via trymail for +dist+/afs/andrew.cmu.edu/usr1/ota/space/space.dl@andrew.cmu.edu (->+dist+/afs/andrew.cmu.edu/usr1/ota/space/space.dl) (->ota+space.digests)
  4.           ID </afs/andrew.cmu.edu/usr1/ota/Mailbox/UZ=kF:y00VcJQHgE5c>;
  5.           Sun,  8 Oct 89 04:24:44 -0400 (EDT)
  6. Message-ID: <kZ=kEe200VcJAHeU44@andrew.cmu.edu>
  7. Reply-To: space+@Andrew.CMU.EDU
  8. From: space-request+@Andrew.CMU.EDU
  9. To: space+@Andrew.CMU.EDU
  10. Date: Sun,  8 Oct 89 04:24:11 -0400 (EDT)
  11. Subject: SPACE Digest V10 #123
  12.  
  13. SPACE Digest                                     Volume 10 : Issue 123
  14.  
  15. Today's Topics:
  16.         C routines for orbital prediction, and a query
  17.         Re: X-30, Space Station Strangles NASP
  18.     Meteorites (was: Re: What to do with the $30 billion)
  19. ----------------------------------------------------------------------
  20.  
  21. Date: 4 Oct 89 12:35:46 GMT
  22. From: ecsvax.uncecs.edu!dukeac!tcamp@mcnc.org  (Ted A. Campbell)
  23. Subject: C routines for orbital prediction, and a query
  24.  
  25.  
  26. [My apologies if this appears twice on your system, but I'm 
  27. reasonably sure that the first time I posted it, it 
  28. didn't make it out.  t.c.]
  29.  
  30. Some readers of sci.astro and sci.space may have seen the Space 
  31. Flight Simulator, version 0.02, which was posted to comp.binaries.ibm.pc.  
  32. If I ever make up my mind just how to do it, I may release 
  33. the source code for the whole thing.  In the meantime, I'm going to 
  34. post the following subroutines to sci.astro and sci.space for 
  35. two reasons (a) someone might actually find them useful, and 
  36. (b) someone might want to advise me on a technical point that 
  37. is apparently beyond the abilities of this bear-of-very-little-brain.  
  38.  
  39. The following code contains three critical subroutines (in C).  
  40. or_init() sets up orbital parameters based on an apoapsis 
  41. and periapsis (in kilometers).  It presumes that the following 
  42. globals have already been set:  
  43.  
  44.     (double) fo_mass    the mass of the orbital focus
  45.                 in grams (earth = 5.98e27)
  46.  
  47.     (double) fo_radius    the radius of the orbital focus 
  48.                 in kilometers
  49.  
  50.     (double) fo_siday    the sidereal period of the orbital 
  51.                 focus in seconds 
  52.  
  53. The include file "sfsm.h" defines a structure for orbits, the elements
  54. of which should be obvious.  sfs_orbits[ n ] is an array of pointers to orbit 
  55. structures, so that sfs_orbits[ 0 ]->period would give the orbital 
  56. period in seconds for orbit 0, etc.  
  57.  
  58. The subroutine or_ssp() calculates the subsatellite point for an
  59. orbit, given the time (in seconds).  It also calculates the 
  60. current altitude.  It, in turn, calls or_kep(), which solves 
  61. Kepler's equation.  These routines may be useful for anyone 
  62. wanting to build C applications with orbital calculations.  
  63. The subroutines are based on pseudocode given in Martin Davidoff's
  64. Satellite Experimenters Handbook (published by ARRL for Radio Amateurs 
  65. using the OSCAR and other Amateur satellites), but the C implementation is
  66. my own.  
  67.  
  68. PROBLEM:  These routines seem to work well enough for orbits with 
  69. inclinations near 0 and near 180 degrees (by the way -- all angles 
  70. in these routines are in radians).  But for orbits with inclinations 
  71. near 90 degrees, there is a grotesque skewing that occurs (at the 
  72. highest/lowest points in the orbit).  As of yet, I haven't been able to figure 
  73. what causes this skewing.  
  74.  
  75. Of course, it's possible that the algorithms are correct, and that 
  76. the skewing occurs in other routines that manipulate their output in 
  77. various ways.  On the other hand, the skewing occurs so frequently 
  78. that I suspect the culprit is here somewhere.
  79.  
  80. email:  tcamp@dukeac.ac.duke.edu
  81.  
  82. /****************************************************************
  83.  
  84.     SFS_OR.C        Orbital calculation routines for SFS
  85.  
  86.             Copyright (c) 1989, Ted A. Campbell
  87.  
  88. ****************************************************************/
  89.  
  90. #include "math.h"
  91. #include "time.h"
  92. #include "vdi.h"
  93. #include "sfsm.h"
  94.  
  95. extern    double or_lfix();
  96.  
  97. /****************************************************************
  98.  
  99.     or_init()        Set up an orbit
  100.  
  101. ****************************************************************/
  102.  
  103. or_init( orbit, ap, aa )
  104.     int orbit;
  105.     double aa, ap;
  106.     {
  107.  
  108.     /***    Calculate semimajor axis        */
  109.  
  110.     sfs_orbits[ orbit ]->semimajor = ( aa + ap + ( 2 * fo_radius )) / 2;
  111.  
  112.     /***    Calculate orbital center - focus distance       */
  113.  
  114.     sfs_orbits[ orbit ]->dist = sfs_orbits[ orbit ]->semimajor - ( fo_radius + ap );
  115.  
  116.     /***    Calculate semiminor axis        */
  117.  
  118.     sfs_orbits[ orbit ]->semiminor = sqrt( (double) ( sfs_orbits[ orbit ]->semimajor * sfs_orbits[ orbit ]->semimajor ) - (double) ( sfs_orbits[ orbit ]->dist * sfs_orbits[ orbit ]->dist ));
  119.  
  120.     /***    Calculate orbital eccentricity  */
  121.  
  122.     sfs_orbits[ orbit ]->eccentricity = sqrt( 1.0 - (( sfs_orbits[ orbit ]->semiminor / (double) sfs_orbits[ orbit ]->semimajor ) * ( sfs_orbits[ orbit ]->semiminor / (double) sfs_orbits[ orbit ]->semimajor )) );
  123.  
  124.     /***    Calculate orbital period        */
  125.  
  126.     sfs_orbits[ orbit ]->period = sqrt ( ( ( 4.0 * ( PI * (double) PI ) ) / ( UGC * fo_mass ) )
  127.            * ( sfs_orbits[ orbit ]->semimajor * sfs_orbits[ orbit ]->semimajor * sfs_orbits[ orbit ]->semimajor ) );
  128.  
  129.     /***    Calculate the increment factor of longitude of 
  130.         the ascending node.  This factor must be multiplied
  131.         by a time factor (orbital period or time into orbit, 
  132.         in seconds.  See Davidoff, pp. 8-9 and 8-10, and 
  133.         formula 8.14.     */
  134.  
  135.     if ( fo_siday == 0.0 )
  136.         {
  137.         sfs_orbits[ orbit ]->lif = 0;
  138.         }
  139.     else
  140.         {
  141.         sfs_orbits[ orbit ]->lif = (2 * PI) / (double) fo_siday;
  142.         }
  143.  
  144.     }
  145.  
  146. /****************************************************************
  147.  
  148.     or_kep()        Solve Kepler's equation
  149.  
  150.     Globals utilized:  
  151.  
  152.     sfs_orbits[ orbit ]->eccentricity Eccentricity of the orbital ellipse
  153.  
  154.     sfs_orbits[ orbit ]->semimajor    Semimajor axis of the orbital ellipse (km)
  155.  
  156.     sfs_orbits[ orbit ]->period       Orbital period (seconds)
  157.  
  158.     Inputs:
  159.  
  160.     t               Time from periapsis in seconds ( 0 < t < sfs_orbits[ orbit ]->period )
  161.  
  162.     Output:
  163.  
  164.     theta           Angle between periapsis, geocenter, and
  165.             current position (theta).  
  166.  
  167.     r               Distance from satellite to focal center,
  168.             in kilometers
  169.  
  170. ****************************************************************/
  171.  
  172. or_kep( orbit, t, theta, r )
  173.     int orbit;
  174.     long t;
  175.     double *theta;
  176.     long *r;
  177.     {
  178.     double z, z3, E;
  179.     z  = 2.0 * PI * ( (double) t / (double) sfs_orbits[ orbit ]->period );
  180.     E  = z;
  181.  
  182.     do
  183.         {
  184.         z3 = ( E - ( sfs_orbits[ orbit ]->eccentricity * sin( E )) - z ) / ( 1 - ( sfs_orbits[ orbit ]->eccentricity * cos( E )));
  185.         E  = E - z3;
  186.         }
  187.     while ( fabs( z3 ) > 0.000000001 );
  188.  
  189.     *theta = PI;
  190.     if ( E != PI )
  191.         {
  192.         *theta = 2.0 * atan( sqrt( ( 1 - sfs_orbits[ orbit ]->eccentricity ) / ( 1 + sfs_orbits[ orbit ]->eccentricity )) * tan( E / 2.0 ));
  193.         if ( E > PI )
  194.             {
  195.             *theta = ( (double) 2.0 * PI ) + *theta;
  196.             }
  197.         }
  198.  
  199.     *r = sfs_orbits[ orbit ]->semimajor * ( 1 - sfs_orbits[ orbit ]->eccentricity * sfs_orbits[ orbit ]->eccentricity ) / ( 1 + sfs_orbits[ orbit ]->eccentricity * cos( *theta ));
  200.     return 1;
  201.     }
  202.  
  203.  
  204. /****************************************************************
  205.  
  206.     or_ssp()    Calculate Subsatellite Point
  207.  
  208.     This function utilizes available data on the satellite 
  209.     to return the subsatellite point, that is, the point 
  210.     on the surface of the orbital focus directly under 
  211.     the satellite at a given moment.  
  212.  
  213. ****************************************************************/
  214.  
  215. or_ssp( orbit, t, latitude, longitude, r, n )
  216.     int orbit;
  217.     long t, *r, *n;
  218.     double *longitude, *latitude;
  219.     {
  220.     static long _r, tp;
  221.     static double theta;
  222.     double n1, n2, n4, E, lan;
  223.     long t1;
  224.  
  225.     *n = t / sfs_orbits[ orbit ]->period;                  /* return number of current orbit */
  226.     t1 = t - ( sfs_orbits[ orbit ]->period * ( *n ) );     /* get seconds into this orbit */
  227.     if ( t1 == 0 )
  228.         {
  229.         t1 = 1;
  230.         }
  231.  
  232.     or_kep( orbit, t1, &theta, &_r );
  233.     *r = _r;
  234.  
  235.     /***    Calculate the longitude of the ascending node at 
  236.         the beginning of this orbit.  See Davidoff, p. 8-9
  237.         and 8-10, and equations 8.13a, 8.13b, and 8.14.    */
  238.  
  239.     lan = sfs_orbits[ orbit ]->lon_an - ( sfs_orbits[ orbit ]->lif * sfs_orbits[ orbit ]->period * *n );
  240.  
  241.     /***    Calculate the latitude of the SSP.  See Davidoff, 
  242.         pp. 8-13 - 8-15, esp. equation 8.20.        */
  243.  
  244.     *latitude = asin( sin( sfs_orbits[ orbit ]->inclination ) * sin( theta + sfs_orbits[ orbit ]->arg_per ));
  245.  
  246.     if ( ( sfs_orbits[ orbit ]->inclination * RAD_DEG >= 0 ) && ( sfs_orbits[ orbit ]->inclination * RAD_DEG < 90 ) )
  247.         {
  248.         n2 = 1;
  249.         }
  250.     else
  251.         {
  252.         n2 = 0;
  253.         }
  254.     if ( ( ( *latitude ) * RAD_DEG ) < 0.0 )
  255.         {
  256.         n4 = 1;
  257.         }
  258.     else
  259.         {
  260.         n4 = 0;
  261.         }
  262.  
  263.     if ( ( sfs_orbits[ orbit ]->arg_per * RAD_DEG > 180 ) && ( sfs_orbits[ orbit ]->arg_per * RAD_DEG <= 540 ) )
  264.         {
  265.         n1 = 1;
  266.         }
  267.     else
  268.         {
  269.         n1 = 0;
  270.         }
  271.  
  272.     E = 2 * atan( (double) pow( ( 1 - sfs_orbits[ orbit ]->eccentricity ) / ( 1 + sfs_orbits[ orbit ]->eccentricity ), (double) 0.5 )
  273.               * tan( sfs_orbits[ orbit ]->arg_per / 2.0 )) + ( 2.0 * PI * n1);
  274.     tp = ( sfs_orbits[ orbit ]->period / ( 2.0 * PI ) ) * ( E - sin( E ) );
  275.  
  276.     *longitude = or_lfix( lan 
  277.              - pow( (double) -1.0, n2 + n4 )
  278.              * acos( cos( theta + sfs_orbits[ orbit ]->arg_per ) / cos( *latitude ) )
  279.              - sfs_orbits[ orbit ]->lif * (double) t1
  280.              - sfs_orbits[ orbit ]->lif * (double) tp );
  281.     }
  282.  
  283. double
  284. or_lfix( l )
  285.     double l;
  286.     {
  287.     double l1;
  288.     l1 = l;
  289.     while ( l1 < ( 0 - PI ) )
  290.         {
  291.         l1 += ( 2.0 * PI );
  292.         }
  293.     while ( l1 > PI )
  294.         {
  295.         l1 -= ( 2.0 * PI );
  296.         }
  297.     return l1;
  298.     }
  299.  
  300. ------------------------------
  301.  
  302. Date: 5 Oct 89 22:23:37 GMT
  303. From: sumax!quick!srg@beaver.cs.washington.edu  (Spencer Garrett)
  304. Subject: Re: X-30, Space Station Strangles NASP
  305.  
  306. In article <SHAFER.89Oct4080229@drynix.dfrf.nasa.gov>, shafer@elxsi.dfrf.nasa.gov (Mary Shafer) writes:
  307. -> 
  308. -> The shuttle comes in from the north to north east.  It comes "feet
  309. -> dry" at about Mach 7 and 145K ft, it's overhead at Edwards at Mach 1
  310. -> at about 40K ft.  It does a HAC (Heading Alignment Circle) to put it
  311. -> on the runway heading (usually 17 or 22), so essentially the pattern
  312. -> is about a 270 teardrop with a longish final.  I think it's a 20 deg
  313. -> glidepath, with a fairly short flair.  Final is flown at 285 KEAS,
  314. -> gear deployed at 275 KEAS, touchdown at 185 KEAS.  (I'm taking these
  315. -> figures from Young and Crippen's 1981 SETP paper on STS-1, so the
  316. -> speeds may not be exact for any given mission, but they're about
  317. -> right.)
  318.  
  319. Um, I believe that should be *towards* N-NE.  Equatorial orbits are
  320. done in the same direction as the Earth's rotation, so they come in
  321. from the W-SW.  (It costs an extra 2000 mph in acceleration to go
  322. the other way.)  Flights employing polar orbits could come in from
  323. either N or S, but I don't think any have yet been made in a Shuttle.
  324. (And may never - I think they need the pad at Vandenburg for polar
  325. launches, and I think that's been officially abandoned.)  Of course,
  326. after that 270 teardrop I believe they do land towards the S-SW,
  327. (ie - from N-NE) so maybe that's what you were talking about.
  328.  
  329. The Shuttle uses *two* glidepaths on final.  They fly most of the
  330. approach on a 17 degree (as I recall) glideslope, then make an
  331. abrupt pitch up to the normal 3 degree glideslope.  It looks
  332. like the two intersect right off the end of the runway (where
  333. "right off" may be a mile or two at these speeds) and they only
  334. spend a few seconds (10 or 15?) on the 3 degree slope before
  335. starting the flare.
  336.  
  337. Now for the questions!  Is "coming feet dry" the same as extending
  338. the landing gear?  I'm pretty sure I remember seeing the gear pop
  339. out *after* the flare, and I can't imagine having to design gear
  340. (much less gear doors) that could handle Mach 7!
  341. And what's the "E" in KEAS?  Surely Edwards doesn't have its own
  342. standard of measurement! :-}
  343.  
  344. ------------------------------
  345.  
  346. Date: Thu, 5 Oct 89 10:20:13 PDT
  347. From: Peter Scott <pjs@aristotle.Jpl.Nasa.Gov>
  348. Subject: Meteorites (was: Re: What to do with the $30 billion)
  349.  
  350. jarvis.csri.toronto.edu!utgpu!utzoo!henry@rutgers.edu  (Henry Spencer) writes:
  351.  
  352. >In article <2153@hydra.gatech.EDU> ccsupos@prism.gatech.EDU (SCHREIBER, O. A.) writes:
  353. >>What is the story about the meteorites findings in Antartica? ...
  354. >
  355. >Well, meteorites fall on Antarctica as they do anywhere else.  Unsurprisingly,
  356. >most of the ones that fall there end up imbedded in ice.  Ice flows, slowly.
  357. >Most of it eventually ends up melting into the ocean, which is not very
  358. >useful.  However, there are a few areas in Antarctica which are "sinks"
  359. >for ice, where flowing ice runs up against rock, is pushed upward, and gets
  360. >steadily eroded by wind.  Anything tough imbedded in that ice ends up sitting
  361. >on the surface.  The result is that meteorites collect on the surface there.
  362.  
  363. Ditto in the Arctic, apparently.  One of the largest meteorites ever recovered
  364. (now on display in the New York Museum of Natural History) was discovered by a
  365. polar explorer wondering where the Eskimo were getting metal from.  They described
  366. a rock that they were chipping it off.  It took two expeditions to finally get
  367. the meteorite (about the size of a VW) and several people-sized companions to the
  368. U.S.  Presumably they compensated the Eskimo with a set of Ginsu steak knives... :-)
  369.  
  370. Peter Scott (pjs@grouch.jpl.nasa.gov)
  371.  
  372. ------------------------------
  373.  
  374. End of SPACE Digest V10 #123
  375. *******************
  376.